PipPenn Data Analysis

Exploration of the Time Series Data

Currently we have data from 10 subjects. Here we only show the examples from Patient 52.

The code for activity types are

##  [1] "Walk"                               "Open Door"                         
##  [3] "Turn Key"                           "Wash Hands"                        
##  [5] "Look around in bag or purse"        "Open Water Bottle and Take a Drink"
##  [7] "Pour Water into a Cup"              "Eat From Bag"                      
##  [9] "1. Patient Choice"                  "2. PatientChoice"                  
## [11] "3. Patient Choice"                  "4. Med Taking"                     
## [13] "5. Med Taking"                      "6. Med Taking"                     
## [15] "7. Med Taking"                      "8. Med Taking"                     
## [17] "9. Med Taking"                      "10. Med Taking"                    
## [19] "11. Med Taking"                     "12. Med Taking"                    
## [21] "13. Med Taking"                     "14. Med Taking"                    
## [23] "15. Med Taking"

Consider “8.Med taking” chunk of activities for analysis.

## look at 8.med taking
med_taking_8 <- subset(dat, Activitylabel == 16)
#dim(med_taking_8)
#min(med_taking_8$time)
med_taking_8$mini_activity <- factor(ifelse(med_taking_8$time > (min(med_taking_8$time)+2)
                                            & med_taking_8$time < (min(med_taking_8$time)+4),
                                            "Open Bottle","Others"), levels = c("Others","Open Bottle"))
## check mini activity
table(med_taking_8$mini_activity)

## check correlation
corrplot(cor(med_taking_8[,c("accel_x_coord","accel_y_coord","accel_z_coord",
                         "rotate_x_coord","rotate_y_coord", "rotate_z_coord")]))

## talk about this with Haochang
plot(med_taking_8$time,med_taking_8$accel_x_coord, type = "l")
lines(rollmean(med_taking_8$accel_x_coord, k = 3), col = "red")
lines(lowess(med_taking_8$time,med_taking_8$accel_x_coord, f = 0.000001), col = "green") 
legend("topleft",                                               # Add legend to line plot
       col = c("red", "green"),
       lty = 1,
       lwd = 2,
       cex = 0.5,
       c("Moving Average, k = 3", "Lowess Smooth, f = 0.000001"))

## fit logistic model
glm_open_bottle <- glm(mini_activity ~ accel_x_coord + accel_y_coord + accel_z_coord +
                         rotate_x_coord + rotate_y_coord + rotate_z_coord + DeviceName,
                         family = "binomial", data = med_taking_8)
summary(glm_open_bottle)
pscl::pR2(glm_open_bottle)

#### perform PCA on 3 accel and 3 rotate variables
predictors <- med_taking_8[,c("accel_x_coord","accel_y_coord","accel_z_coord",
                              "rotate_x_coord","rotate_y_coord", "rotate_z_coord")]
predictors.pca <- prcomp(predictors, center=TRUE, scale.=TRUE)
summary(predictors.pca)
pcs <- predictors.pca$x ## the PCs for regression

## check pc plots
fviz_pca_var(predictors.pca,axes = c(1, 2))
fviz_pca_var(predictors.pca,axes = c(3, 4))
fviz_pca_var(predictors.pca,axes = c(5, 6))

## combine pcs to df
med_taking_8 <- cbind(med_taking_8,pcs)

## perform PCR
glm_pcr <- glm(mini_activity ~ PC1 + PC2 + PC3 + DeviceName, family = "binomial", data = med_taking_8)
summary(glm_pcr)
pscl::pR2(glm_pcr)

#### SVM
svm.df <- med_taking_8[, c("mini_activity","accel_x_coord","accel_y_coord","accel_z_coord",
                           "rotate_x_coord","rotate_y_coord", "rotate_z_coord")]
## find optimal cost of misclassification
# tune.out <- tune(svm, mini_activity ~ .,
#                  data = svm.df,
#                  kernel = "polynomial",
#                  ranges = list(cost = c(0.001, 0.01, 0.1, 1, 5, 10)) 
#                  ) # cost = 1
# ## extract the best model
# (bestmod <- tune.out$best.model)

## Fit Support Vector Machine model to data set
svmfit <- svm(mini_activity ~ ., data = svm.df, kernel = "polynomial", cost = 1)
summary(svmfit)
## Create a table of misclassified observations
ypred <- predict(svmfit, svm.df)
(misclass <- table(predict = ypred, truth = svm.df$mini_activity)) ## NOPE

I also study the movelet package. Thanks to Jiawei example, I was able to put together a workflow and this is what I got:

STEP 1: create an AcceleData using the Acceleration_Create

STEP 2: create chapters of movelets from the AcceleData above using Movelet_Create. Now you have your reference chapters. Next step is prediction of unlabeled chapters

STEP 3: Predict label for each unlabeled chapter of movelets using Movelet_Pred

STEP 4: Validate results using Movelet_Validation

Example:

data.Sample5 <- Acceleration_Create("data_label.csv",subjectName="Sample5",frequency=30)
Ch.sedentary_act.Sample5 <- Movelet_Create(5701,5761,data.Sample5, "sedentary_act")
pred <- Movelet_Pred(data.Sample5,Act.Names=c("sedentary_act"),vote=TRUE)
validation <- Movelet_Validation(data.Sample5, Prediction = pred, Act.Names=c("sedentary_act"))

Try movelet with patient 52

# table(dat$Wrist_Location)
# x = dat$Wrist_Location == "left      "
# test <- ave(x, cumsum(!x), FUN = cumsum)
# table(test)
# which(test > 1)
# x.rm.left <- dat$Wrist_Location[-which(test > 1)]
# table(x.rm.left)

Using 6 PCs for all 15 variables of both hands

## 
##  left right 
## 23300 23300

## extract 3 pcs for the left hand
vars_left <- nodup_dat[nodup_dat$Wrist_Location=="left",
                       c("accel_x_coord","accel_y_coord","accel_z_coord",
                        "rotate_x_coord","rotate_y_coord", "rotate_z_coord",
                        "gravity_x_coord","gravity_z_coord","gravity_y_coord",
                        "magneto_x_coord","magneto_y_coord","magneto_z_coord",
                        "attitude_x_coord","attitude_y_coord","attitude_z_coord")]
varsleft.pca <- prcomp(vars_left, center=TRUE, scale.=TRUE)
#summary(varsleft.pca)
pcs_left <- varsleft.pca$x[,1:3]

## extract 3 pcs for the right hand
vars_right <- nodup_dat[nodup_dat$Wrist_Location=="right",
                       c("accel_x_coord","accel_y_coord","accel_z_coord",
                        "rotate_x_coord","rotate_y_coord", "rotate_z_coord",
                        "gravity_x_coord","gravity_z_coord","gravity_y_coord",
                        "magneto_x_coord","magneto_y_coord","magneto_z_coord",
                        "attitude_x_coord","attitude_y_coord","attitude_z_coord")]
varsright.pca <- prcomp(vars_right, center=TRUE, scale.=TRUE)
#summary(varsright.pca)
pcs_right <- varsright.pca$x[,1:3]

pcs_combined <- cbind(pcs_left,pcs_right)
colnames(pcs_combined) <- c("PC1_left","PC2_left","PC3_left","PC1_right","PC2_right","PC3_right")

# dat$Time <- c(1:nrow(dat))
# table(dat$Wrist_Location)

## combine data with Pcs
ActivityName <- nodup_dat$ActivityName[nodup_dat$Wrist_Location=="left"]
movelet.dat <- data.frame(pcs_combined, ActivityName)
movelet.dat$ActivityName <- gsub(" ", "_", movelet.dat$ActivityName)

## Reorder movelet into left hand then right hand for each activity
#movelet.dat <- movelet.dat[order(movelet.dat$ActivityName,movelet.dat$Wrist_Location),]

## create accel data 
patient.52 <- Acceleration_Create2(movelet.dat,subjectName="patient52",
                                   n_axes = 6,frequency=50) # 50 obs per second
## train pick up and open a pill bottle chapter
## 9._Med_Taking is about 24 seconds, about 1200 observations, so 1 sec = 50 obs. 
from <- min(which(patient.52$Label == "9._Med_Taking")) + 200 ## 12:11:42 picks up pill bottle
to <- min(which(patient.52$Label == "9._Med_Taking")) + 400 ## 112:11:46
Ch.Pick_up_and_open_pill_bottle.patient52 <- Movelet_Create2(FROM = from ,TO = to, 
                                    DATA = patient.52, n_axes = 6,
                                    ChapName = "Pick_up_and_open_pill_bottle")
# ## train pour pill chapter
# from <- min(which(patient.52$Label == "9._Med_Taking")) + 800 ## 12:11:46
# to <- min(which(patient.52$Label == "9._Med_Taking")) + 900 ## 12:11:47
# Ch.Pour_pill.patient52 <- Movelet_Create2(FROM = from ,TO = to,
#                                     DATA = patient.52, n_axes = 6,
#                                     ChapName = "Pour_pill")

## train take pill chapter 
from <- min(which(patient.52$Label == "9._Med_Taking")) + 400 ## 12:11:46 
to <- min(which(patient.52$Label == "9._Med_Taking")) + 500 ## 112:11:48
Ch.Pour_pill_return_bottle_and_take_pill.patient52 <- Movelet_Create2(FROM = from ,TO = to, 
                                    DATA = patient.52, n_axes = 6,
                                    ChapName = "Pour_pill_return_bottle_and_take_pill")

## train Wash_Hands chapter
from <- min(which(patient.52$Label == "Wash_Hands"))
to <- max(which(patient.52$Label == "Wash_Hands"))
Ch.Wash_Hands.patient52 <- Movelet_Create2(FROM = from ,TO = to,
                                    DATA = patient.52, n_axes = 6,
                                    ChapName = "Wash_Hands")
## train Look_around_in_bag_or_purse chapter
from <- min(which(patient.52$Label == "Look_around_in_bag_or_purse"))
to <- max(which(patient.52$Label == "Look_around_in_bag_or_purse"))
Ch.Look_around_in_bag_or_purse.patient52 <- Movelet_Create2(FROM = from ,TO = to,
                                    DATA = patient.52, n_axes = 6,
                                    ChapName = "Look_around_in_bag_or_purse")
## test data
test.dat <- subset(movelet.dat, ActivityName %in% c("Walk","Turn_Key","Open_Door","Eat_From_Bag","6._Med_Taking","7._Med_Taking","4._Med_Taking","5._Med_Taking","8._Med_Taking","10._Med_Taking","11._Med_Taking","12._Med_Taking"))
test.dat.patient.52 <- Acceleration_Create2(test.dat,subjectName="patient52",n_axes = 6,frequency=50) 
train.label=c("Pick_up_and_open_pill_bottle","Pour_pill_return_bottle_and_take_pill","Wash_Hands","Look_around_in_bag_or_purse")

#expand the prediction activity labels to include all trained chapters
pred <- Movelet_Pred2(test.dat.patient.52, Act.Names=train.label, vote=TRUE)
# Save it
#saveRDS(pred, file = "patient52_6pcs_both_separated_wrists_pred.rds")
saveRDS(pred, file = "patient52_8_meds.rds")

# load it
patient52_8_meds <- readRDS(file = "patient52_8_meds.rds")
patient52_8_meds$Pred2 <- recode(patient52_8_meds$Pred, `1`="Pick_up_and_open_pill_bottle", `2`="Pour_pill_return_bottle_and_take_pill", `3`="Wash_Hands", `4`="Look_around_in_bag_or_purse")

## add label for trained activities
test.dat.patient.52$Label2 <- test.dat.patient.52$Label
test.dat.patient.52$Label2[c(
 min(which(test.dat$ActivityName == "4._Med_Taking")):
 min(which(test.dat$ActivityName == "4._Med_Taking")+ 700),
 min(which(test.dat$ActivityName == "5._Med_Taking")):
 min(which(test.dat$ActivityName == "5._Med_Taking")+ 250),
 min(which(test.dat$ActivityName == "6._Med_Taking")):
 min(which(test.dat$ActivityName == "6._Med_Taking")+ 250),
 min(which(test.dat$ActivityName == "7._Med_Taking")):
 min(which(test.dat$ActivityName == "7._Med_Taking")+ 450),
 min(which(test.dat$ActivityName == "8._Med_Taking")):
 min(which(test.dat$ActivityName == "8._Med_Taking")+ 400),
 min(which(test.dat$ActivityName == "10._Med_Taking")):
 min(which(test.dat$ActivityName == "10._Med_Taking")+ 500),
 min(which(test.dat$ActivityName == "11._Med_Taking")):
 min(which(test.dat$ActivityName == "11._Med_Taking")+ 300),
 min(which(test.dat$ActivityName == "12._Med_Taking")):
 min(which(test.dat$ActivityName == "12._Med_Taking")+ 300)
 )] <- "Pick_up_and_open_pill_bottle"

test.dat.patient.52$Label2[c(
 min(which(test.dat$ActivityName == "4._Med_Taking") + 701):
 min(which(test.dat$ActivityName == "4._Med_Taking") + 1150),
 min(which(test.dat$ActivityName == "5._Med_Taking") + 251):
 min(which(test.dat$ActivityName == "5._Med_Taking") + 500),
 min(which(test.dat$ActivityName == "6._Med_Taking") + 251):
 min(which(test.dat$ActivityName == "6._Med_Taking") + 500),
 min(which(test.dat$ActivityName == "7._Med_Taking") + 451):
 min(which(test.dat$ActivityName == "7._Med_Taking") + 650),
 min(which(test.dat$ActivityName == "8._Med_Taking") + 401):
 min(which(test.dat$ActivityName == "8._Med_Taking") + 650),
 min(which(test.dat$ActivityName == "10._Med_Taking") + 501):
 min(which(test.dat$ActivityName == "10._Med_Taking") + 800),
 min(which(test.dat$ActivityName == "11._Med_Taking") + 301):
 min(which(test.dat$ActivityName == "11._Med_Taking") + 650),
 min(which(test.dat$ActivityName == "12._Med_Taking") + 301):
 min(which(test.dat$ActivityName == "12._Med_Taking") + 400)
 )] <- "Pour_pill_return_bottle_and_take_pill"

temp <- data.frame(cbind(test.dat.patient.52$Label2, patient52_8_meds$Pred2, 
                         patient52_8_meds$Pred))
colnames(temp) <- c("Label","Pred","Pred_num")
temp <- temp[test.dat.patient.52$Label2 %in% c("Pick_up_and_open_pill_bottle","Pour_pill_return_bottle_and_take_pill"), ]

## determine number of med taking activities and their location in the data
train <- rle(as.character(temp$Label))
test <- rle(as.numeric(temp$Pred_num))
rle_diff <- diff(test$values)
Pick_up_and_open_pill_bottle_idx <- which(train$values == "Pick_up_and_open_pill_bottle")

## how many % of correct prediction needed to classify a correct prediction?
cut_off <- seq(0.50,0.7,0.05)
pred_res_acc_median <- c()
correct_med_pred <- c()
#true_label <-  c("Med-taking session","Not med-taking","Part of Med-taking, not a full session","Med-taking session","Med-taking session","Med-taking session","Part of Med-taking, not a full session","Med-taking session", "Med-taking session","Med-taking session","Part of Med-taking, not a full session","Part of Med-taking, not a full session","Med-taking session","Part of Med-taking, not a full session")
for (j in 1:length(cut_off)){
  ## determine the cut off for prediction length
  correct_pred_cut_off <- c()
  for(i in 1:length(Pick_up_and_open_pill_bottle_idx)){
    correct_pred_cut_off[i] <- cut_off[j] * sum(train$lengths[Pick_up_and_open_pill_bottle_idx][i] + train$lengths[Pick_up_and_open_pill_bottle_idx + 1][i])
  }
  correct_pred_cut_off_median <- median(correct_pred_cut_off)
  
  ## determine accuracy of the prediction
  pred_res <- c()
  for (i in which(rle_diff == 1)){
    if(test$values[i] == 1 & test$values[i+1] == 2){
      if(test$lengths[i] + test$lengths[i+1] > correct_pred_cut_off_median){
          pred_res <- c(pred_res,"Med-taking session")
        } else {pred_res <- c(pred_res,"Part of Med-taking, not a full session")}
    } else {pred_res <- c(pred_res,"Not med-taking")}
  }
  #pred_res_acc_median[j] <- round(sum(pred_res == true_label)/length(pred_res),2)
  #correct_med_pred[j] <- sum(which(pred_res == "Med-taking session") == which(true_label == "Med-taking session"))
  correct_med_pred[j] <- length(which(pred_res == "Med-taking session"))
}

data.frame(cut_off,correct_med_pred, med_session = 8)

##   cut_off correct_med_pred med_session
## 1    0.50                7           8
## 2    0.55                7           8
## 3    0.60                6           8
## 4    0.65                6           8
## 5    0.70                6           8

## Calculating accuracy
Pick_up_and_open_pill_bottle_acc <- sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle" & as.character(temp$Pred) %in% "Pick_up_and_open_pill_bottle")/sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle")
Pour_pill_return_bottle_and_take_pill_acc <- sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill" & as.character(temp$Pred) %in% "Pour_pill_return_bottle_and_take_pill")/sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill")

## create 2x2 tables
obs_pick_pred_pick <- sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle" & as.character(temp$Pred) %in% "Pick_up_and_open_pill_bottle")
obs_pick_pred_others <- sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle" & !(as.character(temp$Pred) %in% "Pick_up_and_open_pill_bottle"))
obs_pour_pred_others <- sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill" & !(as.character(temp$Pred) %in% "Pour_pill_return_bottle_and_take_pill"))
obs_pour_pred_pour <- sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill" & as.character(temp$Pred) %in% "Pour_pill_return_bottle_and_take_pill")

## row = obs, col = pred
acc_mat <- matrix(c(obs_pick_pred_pick,obs_pick_pred_others,obs_pour_pred_others,obs_pour_pred_pour), byrow = TRUE, nrow = 2)
#acc_mat

The accuracy for predicting “pick up and open pill bottle” is 0.9303 and the accuracy for predicting “pour pill, return bottle and take pill” is 0.5842.

Another way that we determine accuracy is to determine whether the prediction is med-taking related or not. Based on this definition, we construct a confusion matrix below:

#length(test.dat.patient.52$Label)
#length(patient52_8_meds$Pred2)
temp2 <- data.frame(cbind(test.dat.patient.52$Label[1:length(patient52_8_meds$Pred2)], patient52_8_meds$Pred2))
temp2$True_label <- factor(ifelse(temp2$X1 %in% c("Open_Door","Turn_Key","Eat_From_Bag"), "Not Med Taking","Med Taking"))
temp2$Pred_label <- factor(ifelse(temp2$X2 %in% c("Wash_Hands","Look_around_in_bag_or_purse"), "Not Med Taking","Med Taking"))
conf_mat <- confusionMatrix(data = temp2$Pred_label, reference = temp2$True_label)
conf_mat$table

##                 Reference
## Prediction       Med Taking Not Med Taking
##   Med Taking           8628              0
##   Not Med Taking       1331           1154

conf_mat$overall[1]

##  Accuracy 
## 0.8802304

Try movelet with patient 86

#Load Data
# dat_86_88_89<-read.csv('/Users/quycao/Box/PipPenn/Data/Patients-86-88-89.csv',header = TRUE,stringsAsFactors = FALSE)
# dat_86 <- subset(dat_86_88_89, PatientID == 86)
# saveRDS(dat_86, file = "patient_86.rds")
dat <- readRDS(file = "/Users/quycao/Box/PipPenn/Data/patient_86.rds")
dat$Activitylabel<-match(dat$ActivityName,unique(dat$ActivityName))
#convert time series and reorder data by time
dat$timestamp<-substr(dat$activity_timestamp,12,23)
dat$time<-as.numeric(hms(dat$timestamp))
dat <- dat[sort.int(dat$time, index.return = TRUE)$ix, ]
dat$time<-dat$time-dat$time[1]

## Relabel names for wrist location
dat$Wrist_Location <- ifelse(dat$Wrist_Location =="right     ","right","left")

## remove duplicate left and right wrist to create balance data
dat$Wrist_Location2 <- ifelse(dat$Wrist_Location == "left", 0, 1)
keep_idx <- c()
diff_time <- diff(dat$time)
diff_wrist <- diff(dat$Wrist_Location2)

for (i in 1:length(diff_time)){
  if (diff_time[i] <= 0.01 & diff_wrist[i] != 0)
  {
    keep_idx = c(keep_idx, i, i + 1) 
  }
  i = i + 2
}
nodup_dat <- dat[keep_idx, ]
table(nodup_dat$Wrist_Location)

## 
##  left right 
## 19866 19866

## extract 3 pcs for the left hand
vars_left <- nodup_dat[nodup_dat$Wrist_Location=="left",
                       c("accel_x_coord","accel_y_coord","accel_z_coord",
                        "rotate_x_coord","rotate_y_coord", "rotate_z_coord",
                        "gravity_x_coord","gravity_z_coord","gravity_y_coord",
                        "magneto_x_coord","magneto_y_coord","magneto_z_coord",
                        "attitude_x_coord","attitude_y_coord","attitude_z_coord")]
varsleft.pca <- prcomp(vars_left, center=TRUE, scale.=TRUE)
#summary(varsleft.pca)
pcs_left <- varsleft.pca$x[,1:3]

## extract 3 pcs for the right hand
vars_right <- nodup_dat[nodup_dat$Wrist_Location=="right",
                       c("accel_x_coord","accel_y_coord","accel_z_coord",
                        "rotate_x_coord","rotate_y_coord", "rotate_z_coord",
                        "gravity_x_coord","gravity_z_coord","gravity_y_coord",
                        "magneto_x_coord","magneto_y_coord","magneto_z_coord",
                        "attitude_x_coord","attitude_y_coord","attitude_z_coord")]
varsright.pca <- prcomp(vars_right, center=TRUE, scale.=TRUE)
#summary(varsright.pca)
pcs_right <- varsright.pca$x[,1:3]

pcs_combined <- cbind(pcs_left,pcs_right)
colnames(pcs_combined) <- c("PC1_left","PC2_left","PC3_left","PC1_right","PC2_right","PC3_right")

## combine data with Pcs
ActivityName <- nodup_dat$ActivityName[nodup_dat$Wrist_Location=="left"]
movelet.dat <- data.frame(pcs_combined, ActivityName)
movelet.dat$ActivityName <- gsub(" ", "_", movelet.dat$ActivityName)

## Reorder movelet into left hand then right hand for each activity
#movelet.dat <- movelet.dat[order(movelet.dat$ActivityName,movelet.dat$Wrist_Location),]

## create accel data 
patient.86 <- Acceleration_Create2(movelet.dat,subjectName="patient86",
                                   n_axes = 6,frequency=50) # 50 obs per second
## train pick up and open a pill bottle chapter
## 9._Med_Taking is about 18 seconds, about 900 observations, so 1 sec = 50 obs. 
from <- min(which(patient.86$Label == "9._Med_Taking")) + 100 ## 12:40:39 
to <- min(which(patient.86$Label == "9._Med_Taking")) + 300 ## 12:40:43
Ch.Pick_up_and_open_pill_bottle.patient86 <- Movelet_Create2(FROM = from ,TO = to, 
                                    DATA = patient.86, n_axes = 6,
                                    ChapName = "Pick_up_and_open_pill_bottle")
## train pour pill, return bottle and take pill chapter 
from <- min(which(patient.86$Label == "9._Med_Taking")) + 300 ## 12:40:43
to <- min(which(patient.86$Label == "9._Med_Taking")) + 500 ## 12:11:47
Ch.Pour_pill_return_bottle_and_take_pill.patient86 <- Movelet_Create2(FROM = from ,TO = to, 
                                    DATA = patient.86, n_axes = 6,
                                    ChapName = "Pour_pill_return_bottle_and_take_pill")
## train Wash_Hands chapter
from <- min(which(patient.86$Label == "Wash_Hands"))
to <- max(which(patient.86$Label == "Wash_Hands"))
Ch.Wash_Hands.patient86 <- Movelet_Create2(FROM = from ,TO = to,
                                    DATA = patient.86, n_axes = 6,
                                    ChapName = "Wash_Hands")
## train Look_around_in_bag_or_purse chapter
from <- min(which(patient.86$Label == "Look_around_in_bag_or_purse"))
to <- max(which(patient.86$Label == "Look_around_in_bag_or_purse"))
Ch.Look_around_in_bag_or_purse.patient86 <- Movelet_Create2(FROM = from ,TO = to,
                                    DATA = patient.86, n_axes = 6,
                                    ChapName = "Look_around_in_bag_or_purse")
## test data
test.dat <- subset(movelet.dat, ActivityName %in% c("Turn_Key","Open_Door","Eat_From_Bag","6._Med_Taking","7._Med_Taking","8._Med_Taking",
  "10._Med_Taking","11._Med_Taking","12._Med_Taking"))
test.dat.patient.86 <- Acceleration_Create2(test.dat,subjectName="patient86",n_axes = 6,frequency=50) 
train.label=c("Pick_up_and_open_pill_bottle","Pour_pill_return_bottle_and_take_pill","Wash_Hands","Look_around_in_bag_or_purse")

#expand the prediction activity labels to include all trained chapters
pred <- Movelet_Pred2(test.dat.patient.86, Act.Names=train.label, vote=TRUE)
# Save it
#saveRDS(pred, file = "patient86_6pcs_both_separated_wrists_pred.rds")
saveRDS(pred, file = "patient86_6_meds.rds")

# load it
patient86_6_meds <- readRDS(file = "patient86_6_meds.rds")
patient86_6_meds$Pred2 <- recode(patient86_6_meds$Pred, `1`="Pick_up_and_open_pill_bottle", `2`="Pour_pill_return_bottle_and_take_pill", `3`="Wash_Hands", `4`="Look_around_in_bag_or_purse")

## add label for trained activities
test.dat.patient.86$Label2 <- test.dat.patient.86$Label
test.dat.patient.86$Label2[c(
 min(which(test.dat$ActivityName == "6._Med_Taking")+ 100):
 min(which(test.dat$ActivityName == "6._Med_Taking")+ 250),
 min(which(test.dat$ActivityName == "7._Med_Taking")+ 50):
 min(which(test.dat$ActivityName == "7._Med_Taking")+ 300),
 min(which(test.dat$ActivityName == "8._Med_Taking")+ 150):
 min(which(test.dat$ActivityName == "8._Med_Taking")+ 400),
 min(which(test.dat$ActivityName == "10._Med_Taking")+ 100):
 min(which(test.dat$ActivityName == "10._Med_Taking")+ 300),
 min(which(test.dat$ActivityName == "11._Med_Taking")+ 50):
 min(which(test.dat$ActivityName == "11._Med_Taking")+ 250),
 min(which(test.dat$ActivityName == "12._Med_Taking")+ 50):
 min(which(test.dat$ActivityName == "12._Med_Taking")+ 250)
 )] <- "Pick_up_and_open_pill_bottle"

test.dat.patient.86$Label2[c(
 min(which(test.dat$ActivityName == "6._Med_Taking") + 251):
 min(which(test.dat$ActivityName == "6._Med_Taking") + 500),
 min(which(test.dat$ActivityName == "7._Med_Taking") + 301):
 min(which(test.dat$ActivityName == "7._Med_Taking") + 450),
 min(which(test.dat$ActivityName == "8._Med_Taking") + 401):
 min(which(test.dat$ActivityName == "8._Med_Taking") + 500),
 min(which(test.dat$ActivityName == "10._Med_Taking") + 301):
 min(which(test.dat$ActivityName == "10._Med_Taking") + 750),
 min(which(test.dat$ActivityName == "11._Med_Taking") + 251):
 min(which(test.dat$ActivityName == "11._Med_Taking") + 600),
 min(which(test.dat$ActivityName == "12._Med_Taking") + 251):
 min(which(test.dat$ActivityName == "12._Med_Taking") + 400)
 )] <- "Pour_pill_return_bottle_and_take_pill"

temp <- data.frame(cbind(test.dat.patient.86$Label2, patient86_6_meds$Pred2, 
                         patient86_6_meds$Pred))
colnames(temp) <- c("Label","Pred","Pred_num")
temp <- temp[test.dat.patient.86$Label2 %in% c("Pick_up_and_open_pill_bottle","Pour_pill_return_bottle_and_take_pill"), ]

## determine number of med taking activities and their location in the data
train <- rle(as.character(temp$Label))
test <- rle(as.numeric(temp$Pred_num))
rle_diff <- diff(test$values)
Pick_up_and_open_pill_bottle_idx <- which(train$values == "Pick_up_and_open_pill_bottle")

## how many % of correct prediction needed to classify a correct prediction?
cut_off <- seq(0.50,0.7,0.05)
pred_res_acc_median <- c()
correct_med_pred <- c()
#true_label <-  c("Other activities","Med-taking session","Other activities","Med-taking session","Not Med-taking session","Other activities","Med-taking session","Other activities","Med-taking session","Other activities")
for (j in 1:length(cut_off)){
  ## determine the cut off for prediction length
  correct_pred_cut_off <- c()
  for(i in 1:length(Pick_up_and_open_pill_bottle_idx)){
    correct_pred_cut_off[i] <- cut_off[j] * sum(train$lengths[Pick_up_and_open_pill_bottle_idx][i] + train$lengths[Pick_up_and_open_pill_bottle_idx + 1][i])
  }
  correct_pred_cut_off_median <- median(correct_pred_cut_off)
  
  ## determine accuracy of the prediction
  pred_res <- c()
  for (i in which(rle_diff == 1)){
    if(test$values[i] == 1 & test$values[i+1] == 2){
      if(test$lengths[i] + test$lengths[i+1] > correct_pred_cut_off_median){
          pred_res <- c(pred_res,"Med-taking session")
        } else {pred_res <- c(pred_res,"Part of med-taking session")}
    } else {pred_res <- c(pred_res,"Not med-taking session")}
  }
  #pred_res_acc_median[j] <- round(sum(pred_res == true_label)/length(pred_res),2)
  #correct_med_pred[j] <- sum(which(pred_res == "Med-taking session") == which(true_label == "Med-taking session"))
  correct_med_pred[j] <- length(which(pred_res == "Med-taking session"))
}

data.frame(cut_off,correct_med_pred, med_session = 6)

##   cut_off correct_med_pred med_session
## 1    0.50                4           6
## 2    0.55                4           6
## 3    0.60                2           6
## 4    0.65                2           6
## 5    0.70                1           6

## Calculating accuracy
Pick_up_and_open_pill_bottle_acc <- sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle" & as.character(temp$Pred) %in% "Pick_up_and_open_pill_bottle")/sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle")
Pour_pill_return_bottle_and_take_pill_acc <- sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill" & as.character(temp$Pred) %in% "Pour_pill_return_bottle_and_take_pill")/sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill")

## create 2x2 tables
obs_pick_pred_pick <- sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle" & as.character(temp$Pred) %in% "Pick_up_and_open_pill_bottle")
obs_pick_pred_others <- sum(as.character(temp$Label) %in% "Pick_up_and_open_pill_bottle" & !(as.character(temp$Pred) %in% "Pick_up_and_open_pill_bottle"))
obs_pour_pred_others <- sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill" & !(as.character(temp$Pred) %in% "Pour_pill_return_bottle_and_take_pill"))
obs_pour_pred_pour <- sum(as.character(temp$Label) %in% "Pour_pill_return_bottle_and_take_pill" & as.character(temp$Pred) %in% "Pour_pill_return_bottle_and_take_pill")

## row = obs, col = pred
acc_mat <- matrix(c(obs_pick_pred_pick,obs_pick_pred_others,obs_pour_pred_others,obs_pour_pred_pour), byrow = TRUE, nrow = 2)
#acc_mat

For this result, movelet couldn’t detect the first two medication taking sessions. Specifically, the length of prediction for med taking 1 is (72,57) and for med taking 2 is (69,107,192).

test$lengths[1:8]

## [1]  22 137  75  58  64 148 106 135

test$values[1:8]

## [1] 4 3 2 3 2 4 1 3

The accuracy for predicting “pick up and open pill bottle” is 0.6282 and the accuracy for predicting “pour pill, return bottle and take pill” is 0.3731.